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ABSTRACT: The Metci-heuristic approaches can provide a sufficiently good solution to an 
optimization problem, especially with incomplete or imperfect information and with lower 
computational complexity especially for numerical solutions. This paper presents an enhanced PSO 
(Particle Swarm Optimization) technique for solving the same problem. Through the PSO performs well 
but it may require some more iteration to converge or sometimes many repetitions for the complex 
problems. To overcome these problems an enhanced PSO is presented which utilizes the PSO with 
double chaotic maps to perform irregular velocity updates which forces the particles to search greater 
space for best global solution. Finally the comparison between both algorithms is performed for the 
EOQ problem considering deteriorating items, shortages and partially backlogging. The simulation 
results shows that the proposed enhanced PSO converges quickly and found much closer solution then 
PSO. 
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I. INTRODUCTION 

As of late more consideration is paid to the supply chain systems due to increasing expenses, 
globalization of assets, globalization of economy, expanded verities of items and prerequisites of lower reaction 
time. As a result of substantial variability in the impacting conditions, essential suppositions of the EOQ model 
obliges incessant changes as per the inventory. As of late the thoughtfulness regarding an alternate essential 
perspective known as decaying properties of things has been generally examines for the exact demonstrating of 
supply chain systems. 

Since numerous other essential variables may be considered to more practical model arrangement. 
Numerous analysts have officially proposed the various types of EOQ models for diverse things, like no- 
shortage inventory model with constant speed decay, inventory model with Weibull-distributed deteriorating 
rate, inventory model with discount and partial backordering and so on. 

Finally in this paper the enhanced PSO is tested for the EOQ problem considering inventory model 
with deteriorating items is taken [1], which considered that the demand function is exponentially decreasing and 
the backlogging rate is inversely proportional to the waiting time for the next replenishment. The objective of 
the problem is to minimize the total relevant cost by simultaneously optimizing the shortage point and the length 
of cycle. The rest of the paper is organized as the section 2 presents a brief literature review of related articles 
while in section 3 mathematical modeling for the considered EOQ problem is presented. The section 4 describes 
the PSO and the enhanced PSO algorithms followed by the simulation results in section 5. Finally in the section 
6 the conclusion drawn on the basis of simulation results is presented. 



II. LITERATURE REVIEW 

This section presents some of the important literatures related to this work. Kai-Wayne Chuang et al. 
[2] studied pricing strategies in marketing, with objective to find the optimal inventory and pricing strategies for 
maximizing the net present value of total profit over the infinite horizon. The studied two variants of models: 
one without considering shortage, and the other with shortage. Lianxia Zhao [7] studied an inventory model 
with trapezoidal type demand rate and partially backlogging for Weibull-distributed deterioration items and 
derived an optimal inventory replenishment policy. Liang Yuh Ouyang et al. [1] presented an EOQ inventory 
mathematical model for deteriorating items with exponentially decreasing demand. Their model also handles the 
shortages and variable rate partial backordering which dependents on the waiting time for the next 
replenishment. Jonas C.P. Yu [4] developed a deteriorating inventory system with only one supplier and one 
buyer. The system considers the collaboration and trade credit between supplier and buyer. The objective is to 
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maximize the total profit of the whole system when shortage is completely backordered. The literature also 
discuss the negotiation mechanism between supplier and buyer in case of shortages and payment delay. Kuo- 
Lung Hou et al. [10] presents an inventory model for deteriorating items considering the stock-dependent selling 
rate under inflation and time value of money over a finite planning horizon. The model allows shortages and 
partially backlogging at exponential rate. Ching-Fang Lee et al. Al [14], considered system dynamics to propose 
a new order system for integrated inventory model of supply chain for deteriorating items among a supplier, a 
producer and a buyer. The system covers the dynamics of complex relation due to time evolution. Michal 
Pluhacek et al [15] compared the performance of two popular evolutionary computational techniques 
(particle swarm optimization and differential evolution) is compared in the task of batch reactor geometry 
optimization. Both algorithms are enhanced with chaotic pseudo-random number generator (CPRNG) based on 
Lozi chaotic map. The application of Chaos Embedded Particle Swarm Optimization for PID Parameter Tuning 
is presented in [16], Magnus Erik et al [17] gives a list of good choices of parameters for various optimization 
scenarios which should help the practitioner achieve better results with little effort. 



III. MATHEMATICAL MODELING OF SYSTEM AND PROBLEM FORMULATION 

The mathematical model in this paper is rendered from reference [1] with following notation and 
assumptions. 

Notation: 

C| : Holding cost, ($/per unit)/per unit time. 
c 2 : Cost of the inventory item, $/per unit. 
c 3 : Ordering cost of inventory, $/per order. 
c 4 : Shortage cost, ($/per unit)/per unit time. 
c 5 : Opportunity cost due to lost sales, $/per unit, 
t i : Time at which shortages start. 

T : Length of each ordering cycle. 

W : The maximum inventory level for each ordering cycle. 

S : The maximum amount of demand backlogged for each ordering cycle. 

Q : The order quantity for each ordering cycle. 

Inv(t ) : The inventory level at time t. 

Assumptions: 

1 . The inventory system involves only one item and the planning horizon is infinite. 

2. The replenishment occurs instantaneously at an infinite rate. 

3. The deteriorating rate, 9 (0 < 9 < 1), is constant and there is no replacement or repair of deteriorated units 
during the period under consideration. 

4. The demand rate P(t), is known and decreases exponentially. 

{Ae~^,l{t) > 0 

D , I{t) < 0 

Where A (> 0) is initial demand and A (0 < A < 9) is a constant governing the decreasing rate of the demand. 

5. During the shortage period, the backlogging rate is variable and is dependent on the length of the waiting time 
for the next replenishment. The longer the waiting time is, the smaller the backlogging rate would be. Hence, the 
proportion of customers who would like to accept backlogging at time t is decreasing with the waiting time 
( T — t ) waiting for the next replenishment. To take care of this situation we have defined the backlogging rate 



m = ^ 



to be ■ 



l+scr-t) 



when inventory is negative. The backlogging parameter 8 is a positive constant fq < t < T . 



3.1 Model Formulation 

Here, the replenishment policy of a deteriorating item with partial backlogging is considered. The 
objective of the inventory problem is to determine the optimal order quantity and the length of ordering cycle so 
as to keep the total relevant cost as low as possible. The behavior of inventory system at any time is depicted in 
Figure 1. 
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Replenishment is made at time t = 0 and the inventory level is at its maximum W. Due to both the 
market demand and deterioration of the item, the inventory level decreases during the period [0, t x ] , and 
ultimately falls to zero at t = t l . Thereafter, shortages are allowed to occur during the time interval [t 1 ,T] and 
all of the demand during the period [t 1( T] is partially backlogged. 

As described above, the inventory level decreases owing to demand rate as well as deterioration during 
inventory interval [0, tj. Hence, the differential equation representing the inventory status is given by 

dlnv(t ) 

1- 6Inv(t ) = — Ae At ,0 < t < ti (3-1) 

dt 

with the boundary condition / (0) = W. The solution of equation ( 1) is 

jA.6 ^ 

Inv(t ) = — [e^-Dbi-D — l], 0 < t < tq (3-2) 

So the maximum inventory level for each cycle can be obtained as 

W = Inv( 0) = l e ?~ X)t - !] (3-3) 

During the shortage interval [t 1 ,T], the demand at time t is partly backlogged at the fraction — Thus, the 



differential equation governing the amount of demand backlogged is as below. 

dlnv(t ) D 

— = (3 - 4) 

with the boundary condition I(tf) = 0 . The solution of equation (3.4) can be given by 

D 

Inv(t ) = — {ln[l + 8(T — £)] — InL'pL + S(T — tq)]}, < t <T 



....(3.5) 

Let t = T in (3.5), we obtain the maximum amount of demand backlogged per cycle as follows: 

D 

S = - Inv(T ) = — ln[l + S(T - tq)] (3.6) 



Hence, the order quantity per cycle is given by 

^4 

Q = W + S = \ P IB-X)t i _ 



6 -A 

The inventory holding cost per cycle is 

rtl 



[ e (9-A)ti _ i] + £ ln [! + S(T - tq)] (3.7) 

A 



HC = f c 1 /nr , (t)dt = 
J n 



c 1 A 



6(6 - A) 



3 “'It 1 



e etl - 1 --r(e itl - l) 

A 



.(3.8) 



The deterioration cost per cycle is 
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DC = c 2 [W - 



[ R(t)dt] 
Jn 



— c iA 

The shortage cost per cycle is 



= c 2 [W- f Ae 
J o 

[jZJ 0 (9 " A)tl - 0 - J (1 - e_lt1 )} (3-9) 



SC = c 4 



— I I(t)dt 

'ti 

The opportunity cost due to lost sales per cycle is 
rTr 1 



j I(t)dt = c 4 D { T —^ - ^ln[l + S(T - q)] j ... . (3.10) 



BC 



- 

Jt.A 



1 



Ddt = c 5 D | (T - - ^ln[l + 5(T - q)]J ... (3.11) 



1 + 5(T - t) 

Therefore, the average total cost per unit time per cycle is 

TVC = TVC{t v T) 

= (holding cost + deterioration cost + ordering cost + shortage cost + opportunity cost due to lost sales)/ length 
of ordering cycle 



TVC = - 



c 4 A 



,-At, a et 



T [6(9 - X) 



[ e 8tl - 1 - - (e Atl - l)| + c 2 A 
c 3 d(^ + c s ) T-t 



^(0— A)ti i g— At 1 



e-x 

ln[l + 8(T - tO] 



1 (A(c-i + dc 2 ) r ra n 1 

TVC = - ' * - ( B-X) tl - 1 ] - 



A(c 4 + 9c 2 ) 



T [ B{B -X) 

D(c 4 + Scf) 



+ ■ 



T-U- 



9X 

ln[l + S(T - t 4 )] 



[l — Xt 4 — e /ltl ] + c 3 



.(3.12) 



The objective of the model is to determine the optimal values of t, and T in order to minimize the average total 
cost per unit time, TVC. The optimal solutions tj 1 and T* need to satisfy the following equations: 



dTVC _ 1 f/Kq +0c 2 ) r 
dt 4 -f{ 9 ^ 



(0— X)t\ p—kt i] 



D(c 4 + Sc 5 ) 



1 - 



= 0. 



and 



dTVC _ 1 (D(c 4 + bc 5 ) 



dT 



T 2 [ 



(T — ti)(5ti — 1) 1 

i_^_i + _ ln[1 + s(7 -_ tl) , 



1 + 5(T - t 4 ) 

^(c 4 + 9c 2 ) 



.(3.13) 



0(0 -X) 






c 3 



= 0 ... (3.14) 



IV. PARTICLE SWARM OPTIMIZATION (PSO) 

The PSO algorithm is inspired by the natural swarm behavior of birds and fish. It was introduced by 
Eberhart and Kennedy in 1995 as an alternative to other ECTs, such as Ant Colony Optimization, Genetic 
Algorithms (GA) or Differential Evolution (DE). Each particle in the population represents a possible solution 
of the optimization problem, which is defined by its cost function. In each iteration, a new location (combination 
of cost function parameters) of the particle is calculated based on its previous location and velocity vector 
(velocity vector contains particle velocity for each dimension of the problem). 

The PSO algorithm works by simultaneously maintaining several candidate solutions in the search 
space. During each iteration of the algorithm, each candidate solution is evaluated by the objective function 
being optimized, determining the fitness of that solution. Each candidate solution can be thought of as a particle 
“flying” through the fitness landscape finding the maximum or minimum of the objective function. 

Initially, the PSO algorithm chooses candidate solutions randomly within the search space. It should be 
noted that the PSO algorithm has no knowledge of the underlying objective function, and thus has no way of 
knowing if any of the candidate solutions are near to or far away from a local or global maximum. The PSO 
algorithm simply uses the objective function to evaluate its candidate solutions, and operates upon the resultant 
fitness values. 

Each particle maintains its position, composed of the candidate solution and its evaluated fitness and its 
velocity. Additionally, it remembers the best fitness value it has achieved thus far during the operation of the 
algorithm, referred to as the individual best fitness and the candidate solution that achieved this fitness, referred 
to as the individual best position or individual best candidate solution. Finally, the PSO algorithm maintains the 
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best fitness value achieved among all particles in the swarm, called the global best fitness, and the candidate 
solution that achieved this fitness, called the global best position or global best candidate solution. 

The PSO algorithm consists of just three steps, which are repeated until some stopping condition is 

met: 

1 . Evaluate the fitness of each particle 

2. Update individual and global best fitness’s and positions 

3. Update velocity and position of each particle 

4. Repeat the whole process till the 

The first two steps are fairly trivial. Fitness evaluation is conducted by supplying the candidate solution 
to the objective function. Individual and global best fitness’s and positions are updated by comparing the newly 
evaluated finesses against the previous individual and global best fitness’s, and replacing the best fitness’s and 
positions as necessary. 

The velocity and position update step is responsible for the optimization ability of the PSO algorithm. 
The velocity of each particle in the swarm is updated using the following equation: 

v(i + 1) = w * v(i) + Cx * ( pBest — x(t)) + c 2 * ( gBest — x(t)) (4.1) 

Where: 

v(i + 1) — New velocity of a particle. 

v(i) — Current velocity of a particle. 

c i > c 2 ~ Priority factors. 

pBest — Best solution found by a particle. 

gBest — Best solution found in a population. 

x(t) — Current position of a particle. 

The new position of a particle is then given by (4.2), where x(t + 1) is the new position: 

x(t + 1) = x(t) + v(i + 1) (4.2) 

Each of the three terms ( w * v(i), q * ( pBest — x(t)) and c 2 * ( gBest — x(t)) of the velocity update 
equation have different roles in the PSO algorithm. The first term w is the inertia component, responsible for 
keeping the particle moving in the same direction it was originally heading. The value of the inertial coefficient 
w is typically between 0.8 and 1.2, which can either dampen the particle’s inertia or accelerate the particle in its 
original direction. Generally, lower values of the inertial coefficient speed up the convergence of the swarm to 
optima, and higher values of the inertial coefficient encourage exploration of the entire search space. 

The second term c 1 * ( pBest — x(t)) called the cognitive component, acts as the particle’s memory, 
causing it to tend to return to the regions of the search space in which it has experienced high individual fitness. 

The cognitive coefficient c 1 is usually close to 2, and affects the size of the step the particle takes 
toward its individual best candidate solution pBest. 

The third termc 2 * ( gBest — x(i)), called the social component, causes the particle to move to the 
best region the swarm has found so far. The social coefficient c 2 is typically close to 2, and represents the size of 
the step the particle takes toward the global best candidate solution gBest the swarm has found up until that 
point. 

4.1 Enhanced PSO 

In the Enhanced PSO direct output iterations of the chaotic map were used for the generation of real 
numbers for the main PSO formula that determines new velocity, thus the position of each particle in the next 
iteration (See Rand in equation (4.1) in section 5). The concept of embedding chaotic dynamics into 
evolutionary algorithms is presented in [15] [16]. The two different maps which are used here are 

4.1.1 Logistic map 

The logistic map is a polynomial mapping (equivalently, recurrence relation) of degree 2, shows the 
complex, chaotic behavior from very simple non-linear dynamical equations. Mathematically, the logistic map 
is written 

X n+1 =pX n a~X n ) (4.3) 



4.1.2. lozi map 

The Lozi map is a simple discrete two-dimensional chaotic map. The map equations are given in (17). 

X n+1 = l-x\ X n \ + bY n (4.4a) 

Yn+i =X n (4.4 b) 
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V. SIMULATION AND RESULTS 

To test the PSO and Enhanced PSO algorithms for finding the optimal solution of the problem as 
discussed in section (3), the equation (3.13) and (3.14) is taken as the fitness function, while the values of t 1 and 
T are to be searched for satisfying the fitness function. The other variables are taken as follows: 



Table I 



Variable Name 


Variable Value 


A 


12 


6 


0.08 


8 


2 


A 


0.03 


0 


0.5 


c 2 


1.5 




10 


c 4 


2.5 


c 5 


2 


D 


8 



The exact solution found analytically in [1] was t 1 = 1.4775 and T = 1.8536. Which provides the optimal 
maximum inventory level W = 18.401 units, the optimal order quantity Q = 20.1183 units and the minimum 
average total cost per unit time TV C — $11.1625. 

While the PSO based technique with following configuration- 



Table II 



Variable Name 


Variable Value 


c 4 (PSO Const.) 


2 


c 2 (PSO Const.) 


1 


Wmax 


0.9 


w min 


0.1 


Total Particles 


100 


Maximum Iterations 


100 


p ( Logistic Map ) 


4 (Only For E-PSO) 


k ( Logistic Map ) 


0.63 (Only For E-PSO) 


a ( Lozi Map ) 


1.7 (Only For E-PSO) 


b ( Lozi Map ) 


0.5 (Only For E-PSO) 



The equation (3.12) can be directly selected as objective function for the PSO algorithm hence it avoids the 
complex differentiations required for analytical methods. 

Using the equation (3. 12) as objective function the results obtained by all three PSO as listed below- 



Table III: 



Technique 


h 


T 


TVC 


PSO 


1.4785 


1.8550 


11.1625 


PSOl 


1.4772 


1.8532 


11.1625 


PS 02 


1.4771 


1.8531 


11.1625 


Analytical 


1.4775 


1.8536 


11.1625 



Every algorithm Provides the minimum value of TVC — $11.1625, although at different t, and T. 



I UMER I ISSN: 2249-6645 I 



www.ijmer.com 



I Vol. 5 I Iss.ll Jan. 2015 I 341 



An Enhance PSO Based Approach for Solving Economical Ordered Quantity (EOQ) Problem 




Figure 2: Comparison for the value of objective function (fitness value) at each iteration of all three PSO (PSO 
= standard PSO, PSOl = PSO with Single Chaotic Map, PS02 = PSO with Double Chaotic Map. 



Also if we try to use the objective function by combining the equation (13) and (14) 

objective function = abs(eq( 13)) + abs(eq( 14)) 

Using the above given equation as objective function the results obtained by all three PSO as listed below 



Technique 


U 


T 


obj 


PSO 


1.4776 


1.8538 


0.3417e-3 


PSOl 


1.4776 


1.8537 


0.4098e-3 


PS02 


1.4775 


1.8537 


0. 1943e-3 


Analytical 


1.4775 


1.8536 


0 




0 10 20 30 40 50 60 70 

Iterations 



90 100 



Figure 3: Comparison for the value of objective function (fitness value) at each iteration of all three PSO. The 
best fitness value achieved in the process is 1.943 * 10 -04 by PS02 which almost equal to 0. 



VI. CONCLUSIONS 

This paper presented an enhanced PSO algorithm for efficient solution of classical economic order 
quantity (EOQ) model problem. The simulation results shows that the enhanced PSO not only converges 
quickly about half of the iterations of standard PSO but also provides much closer solution. Hence the presented 
approach can solve the EOQ problem much efficiently than standard PSO and can be preferred for complex 
EOQ models where analytical solution required extensive analysis or normal PSO fails to produce required 
closeness to solution. 
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